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Structure and dynamics of penta-hepta defects in hexago- 
nal patterns is studied in the framework of coupled amplitude 
equations for underlying plane waves. Analytical solution for 
phase field of moving PHD is found in the far field, which gen- 
eralizes the static solution due to Pismen and Nepomnyashchy 
(1993). The mobility tensor of PHD is calculated using com- 
bined analytical and numerical approach. The results for the 
velocity of PHD climbing in slightly non-optimal hexagonal 
patterns are compared with numerical simulations of ampli- 
tude equations. Interaction of penta-hepta defects in optimal 
hexagonal patterns is considered. 

Hexagonal patterns in nonequilibrium extended sys- 
tems are formed as a result of superposition of three plane 
waves oriented at 120° with respect to each other. They 
appear naturally in a large variety of pattern-forming 
systems in fluid dynamics [1], optics [2], chemical kinet- 
ics [3], etc. The most generic defect in hexagonal pattern 
is so called penta-hepta defect (PHD) which is a bound 
state of two dislocations of opposite winding numbers 
on two different wave systems [1,4]. In the paper [5] 
we have demonstrated that the mechanism which pro- 
vides binding force is the synchronization of the phase 
fields in the bulk of the system so the resonant condition 
for phases can be fulfilled everywhere except the core 
of the defect. If the wavenumbers of waves composing 
hexagonal pattern are equal to the onset value of cellular 
instability, PHD stays put. Meanwhile, observations of 
penta-hepta defects in non-Boussinesq Rayleigh-Benard 
convection and in other systems demonstrate that PHDs 
in fact move slowly and eventually annihilate or disap- 
pear at the boundaries [6]. It is conceivable that the mo- 
tion of PHDs is caused by ambient strain due to deviation 
of wavenumbers from onset value or by defect interaction 
with each other. An example of non-optimal hexagonal 
pattern, non-eqilateral (rhombic) hexagons, was consid- 
ered recently in a number of papers [7-9]. As we shall 
see, even equilateral hexagonal patterns with wavenum- 
ber different from the optimal value, produce a driving 
Peach-Kohler force for penta-hepta defect. These non- 
optimal equilateral pattern were considered in [10-12]. 

In this Letter, the motion of PHDs in non-optimal 
hexagonal patterns is considered. Analysis is carried in 
the framework of three coupled equations for complex 
amplitudes of individual waves. These equations con- 
stitute the simplest heuristic model for hexagonal pat- 
terns and play essentially the same role here as com- 
plex Ginzburg-Landau equation for oscillatory media 
or Newell- Whitehead-Segel equation for roll patterns in 
isotropic media. 



The order parameter of the hexagonal pattern is writ- 
ten in the form a = e 1 / 2 X2j=i Bj exp(i(kj + K/)r) + c.c. 
where e is a small parameter characterizing the distance 
above onset, all three |k;| = ko , an onset wavenumber of 
the symmetry-breaking instability, ki points at the x di- 
rection (polar angle <f> = 0), and k2 ) 3 point at <f> = ±27r/3 
respectively. Kj are rescaled corrections to the optimal 
wavevectors k j , Yl^=i^-j = from the resonance con- 
dition. Complex amplitudes Bj are slow functions of 
R = e 1/,2 r and T = et. Bj satisfy the following triplet of 
equations [11]: 

drBj = 0* - Kf)Bj + B*_,B* +1 - (\B 3 \ 2 + 7 |B J -_ 1 | 2 

+ 7 |5 J+1 \ 2 )Bi + (uj • VfB 3 + 2/ A in • V)B 3 . (1) 

Here rij is the unit vector oriented along the wavevector 

of wave j, m = (1, 0), n 2 = (-±, ^), n 3 = {-\, 
index j is defined modulo 3, the coefficient of quadratic 
nonlinearity is rescaled to unity, \i is the rescaled super- 
criticality parameter, 7 is the ratio of the coefficient of 
cubic interaction of rolls of different orientation to the co- 
efficient of cubic self-interaction, and Kj = (rij-Kj). Spa- 
tial gradients are calculated with respect to slow variable 
R, and asterisks denote complex conjugate. The applica- 
bility of these equations to description of real hexagonal 
patterns has been discussed earlier [5,11,13]. 

Eqs.(l) have a family of uniform stationary solutions, 
Bj , of which only equilateral one (Kj = K) can be ex- 
pressed in a simple analytic form 

R o T3 - 1 + V 1 + 4 0» - A ' 2 )(! + 2 T) (9) 

In a general case of different Kj the amplitudes Bj are 
different, and the solutions are known as rhombic pat- 
terns. Stability properties of nonequilateral rhombic pat- 
terns were studied recently in Ref. [7,8] . 

Spontaneously formed hexagonal patterns are usually 
defect-ridden. Various defects have been described in 
the literature (see, for example, [14]). Most of them are 
not stable and either disappear quickly or transform into 
basic penta-hepta defects. As it was mentioned before, 
PHD is a bound state in which two of three modes have 
dislocations with opposite winding numbers. Without 
loss of generality we will consider a particular form of 
penta-hepta defect with positive dislocation in mode 1 
and negative dislocation in mode 3. Mode 1 contains no 
dislocations. Corresponding solution to (1) can be writ- 
ten in the form Aj = Fj(R, </>) e lS ^ R '^\ where R and 



1 



<f> are polar coordinates, § c V9\ds = 0, § c V#2,3^s = 
±27r, 7*12(0) = 0, 7*i_ 3 (oo) = Bo, and C is a closed 
contour encircling the origin. This solution cannot be 
expressed in a closed analytic form. However, in the far 
field where all the amplitudes approach asymptotic value 
Bo, the following solution for the phase fields 6j depend- 
ing only on the polar angle <f> has been found [13]: 

9 1 = (1 - cos 2«/>)y , #2,3 = ±<f>-[\- cos(2^ T y )]y 

(3) 

In order to find the equations of motion for the penta- 
hepta defect in the non-optimal hexagonal pattern, we as- 
sume that PHD moves with a constant velocity V. Trans- 
forming into a moving frame R' = R — VT then yields 
the set of stationary equations for Bj(R') which coincides 
with (1) with 8t replaced by — VVBj. Then we project 
these equations onto its two orthogonal translational 
modes, and {d*Bj} (we choose a new coordinate 

frame (£,77) where £ = X cos ip + Ysint/> is the coordi- 
nate along the defect motion, and rj = Y sin ip — X cos ip 
orthogonal to that, ip is the angle between the direction 
of defect motion and X-axis): 

where I is a mobility tensor of PHD, 

3 3 
% = <El3e5 J f), I vv =(£\d v Btf), 
i=i i=i 

1 3 

i=i 

Ti = 27r[|5 2 °| 2 A' 2 sm(^ - y ) - |B 3 ° | 2 A' 3 sinty + y)], 
T 2 = 2tt[\B° 2 \ 2 K 2 cos(^ - y ) - |5 3 | 2 A' 3 costy + y )]. 

and (...) = Jf...dXdY. All other terms from the r.h.s. 
of eq.(l) vanish upon integration under usual boundary 
conditions at infinity. Here we used the formula 

(d x Bjd Y B*) - c.c. = 2m6j |5°| 2 , (6) 

where 6j is a winding number of the dislocation at the 
particular mode, 61 = 0, and 82,3 = ±1, respectively. 

Equations of motion (4) indicate that the penta-hepta 
defect is driven by superposition of two Peach-Kohler 
forces corresponding to a strain at singular modes 2,3. 
If K'2,,3 = 0, defect does not move. Strain at non-singular 
mode 1 does not enter explicitly in the equations (4), 
however implicitly K\ affects the amplitudes Bj . 

The well known difficulty in treating equations of mo- 
tions for topological defects is that integrals entering 



their mobilities diverge at large distances when station- 
ary solutions are used in the integrands. In fact, with 
the static phase approximation solution (3) components 
of I diverge logarithmically at both small and large R. 
At small R the phase approximation is not valid, and the 
stationary solution of full amplitude equations (1) should 
be employed. More serious problem arizes at large R. Ev- 
idently, one could introduce an ad hoc large-scale cut-off 
due to finite-size effects, and therefore the mobility will 
be logarithmically dependent on the size of the box. This 
may be relevant only for small systems (V Rbox "C 1). 
Another possibility widely considered in the literature 
[15,16] for regular dislocations is to use solutions corre- 
sponding to moving defects [17]. In this case the integrals 
converge and therefore a finite velocity of dislocations can 
be found even in the large box limit (V Rbox ^1)- In the 
spirit of the calculation scheme [16] we assume that all 
Kj are small, therefore velocity of PHD is also small, and 
the solution describing moving defect differs from the sta- 
tionary one only at large distances R ~ V~ x ^> 1, where 
the phase approximation is well justified. Therefore, the 
region of integration for the components of the mobility 
tensor can be split into two parts; in the inner region 
(R < R , where 1 < R < | V| _1 ) the stationary PHD 
solution can be used, and in the outer region (R > Ro) 
the phase approximation can be used to simplify the task 
of finding the moving PHD solution. At R ~ Ro static 
phase approximation solution (3) is applicable, so both 
parts should depend on Ro logarithmically. After adding 
them together radius Ro should drop out. 

The details of this rather involved calculations will be 
published elsewhere [19]. The final expressions for the 
mobility components are relatively simple: 

7 ee = -5o[y H w i y ) - 7rcos2i/>m(w2V)] 

1^ = 5 2 7rsin2t/>ln(w 3 V0, (7) 

where w\ = 1.24, u>2 = 1.14, W3 = 2.00 (we do not need 

I'll)- 

Now we can substitute (5) and (7) into (4). Since we 
already assumed all three Kj small, without further loss 
of accuracy we can take in (5) all three \Bj\ = Bo af- 
ter which it is cancelled. Thus we obtain a set of two 
nonlinear algebraic equations for V and t/>, 

^Vln(wiV) + V\n(w 2 V) cos 2t/> = 

2K 2 sm(-iP + y) + 2/i' 3 sin(^ + y), (8) 

27T 27T 

Vln(w 3 V)sm2il> = 2K 2 cos(-i/> + — ) - 2/i' 3 cos(t/> + — ) 

o o 

In these equations there are only three 0(1) constants 
which (for small Kj) are functions of parameters /i,y 
only. 

Equations (4) have been written for a particular penta- 
hepta defect with positive dislocation in the second 
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mode, negative dislocation in the third mod and no dis- 
location in the first mode. We label it as (0,1,-1). 
Totally, there exist six distinct penta-hepta defects, 
(0, 1, —1), (1, 0, —1), (1, —1, 0), and their mirror images 
(conjugate defects) (0, -1, 1), (-1, 0, 1), (-1, 1, 0). The 
equations for (1, 0, —1), (1, —1, 0) PHDs can be obtained 
from Eqs.(8) by cyclic relabeling of ^1,2,3 an d replacing 
ip by ip ± ^. For a conjugate PHD, the mobility tensor 
remains the same, but the r.h.s. of equations of motion 
change sign. Obviously, for conjugate defects, V* = V 
and ip* = ip + 7r. 

Numerical simulations of the amplitude equations (1) 
were performed using a split-step method. Linear parts 
were integrated using FFT, and nonlinear parts were 
calculated using Euler integrator. Typically we used 
256 2 spatial harmonics with periodic boundary condi- 
tions; physical system size was 100, and the time step was 
chosen 0.1. In all examples described below = 1, 7 = 2. 
As initial condition we take (0,1,-1) defect placed in the 
middle of the integration domain. To diminish an effect 
of periodic boundary conditions we introduce a circular 
ramp at R > 0AL. Detailed results of numerical calcu- 
lations will be given in [19]. Here we present only one 
pertinent example. 

In the Figure 1 the magnitude and angle of the ve- 
locity vector found from (4) are plotted together with 
the results of direct numerical simulations versus K3 for 
K\ = A'2 =0.1. Quantitative comparison indicates that 
both the direction of motion and the magnitude of veloc- 
ity are in a good agreement with the theoretical analysis. 
We checked the importance of finite-size effects by com- 
puting velocity for a smaller system size L = 40. The 
results remained practically the same. 
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FIG. 1. Velocity vector of a PHD (0, 1, —1) as a function of 
wavenumber correction if 2 at K\ = 0.1, K3 =0.1; main panel 

- angle with respect to X-axis, parameters fi = 1, 7 = 2, 
inset - magnitude. Solid lines - theory, Eqs.(8), open circles 

- numerical simulations of (1) with system size L = 40, solid 
circles - same for L = 100 

The theory predicts (and numerical calculations con- 
firm) that the direction of PHD motion strongly depends 
on the combination of the wavenumber corrections at the 



singular modes. On the contrary, wavenumber detuning 
at the non-singular mode K\ only weakly affects V via 
amplitudes B® . Notably, even if all three wavevectors 
are equal (but non-optimal), PHD still moves along the 
X-axis (wavevector of non-singular mode). 

It is tempting to apply the equations of motion (8) 
directly to the interaction of two PHDs. Indeed, when 
two PHDs are far enough, they interact entirely through 
phase perturbations. Each defect distorts the phase field 
and therefore creates slightly non-optimal hexagonal pat- 
tern at the location of another defect. However, if we 
assume that the field is static, and estimate convergence 
velocity for a distance R, we obtain V R ~ 2 which clearly 
is inconsistent with this assumption. Strictly speaking, 
phase field "remembers" whole previous path of the de- 
fects. 

Here we present some of the results of numerical simu- 
lations of interacting penta-hepta defects in an ideal (all 
Kj = 0) pattern. The result of interaction (attraction or 
repulsion) of two PHDs, (8\ , 8\, 63) and (8f , <5|, <$f), de- 
pends on the number N = Y^j=i ^j^j ■ ^ can on ly take 
values of —2, —1, 1, and 2. If N < defects attract each 
other, and in all other cases they repel each other. In the 
Figures 2a, b two families of PHD trajectories are shown 
for several initial positions of defects. Figure 2a illus- 
trates attraction of two conjugate defects, (0, 1, —1) and 
(0, —1, 1) (in this case N = —2 and defects are attracted). 
After collision pairs of dislocations at modes 2 and 3 an- 
nihilate, and thus perfect hexagonal pattern establishes. 
Meanwhile, two PHDs of the same type, (0, 1, —1) (here 
N = 2), repel each other (see Fig. 2b). For two different 
PHDs, (0, 1, -1) and (-1, 0, 1), N = -1, and defects are 
attracted again, however complete annihilation does not 
occur. Instead, conjugate dislocations at mode 3 anni- 
hilate, and remaining dislocations at modes 1 (from the 
first PHD) and 2 (from the second PHD) immediately 
form a new penta-hepta defect, (—1,1,0). As this de- 
fect is alone, and the ambient strain is absent, the defect 
stays put. 
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FIG. 2. Trajectories of interacting penta-hepta defects. 
Open circles indicate positions of the defect cores with time 
interval AT = 2.5. Arrows point toward directions of motion; 
a - (0,1,-1) (solid) and (0,-1,1) (dashed); b - (0,1,-1) 
and (0, 1, —1). For each case, four sets of initial conditions 
are taken (they are labeled 1-4 in the figures). 
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FIG. 3. Distances between interacting penta-hepta defects 
(0, 1, —1) and (0, —1, 1) versus time, a - linear coordinates, b 
- logarithmic coordinates. Labels 1-4 correspond to different 
initial positions of the defects, as shown in Fig. 2a. 

In the Figure 3a the distance between two cores is 
shown as a function of time for interacting (0,1,-1) 
and (0,-1,1) defects. The rate of convergence varies 
with initial positions of defects. Attraction is strongest 
for defects aligned along Y-axis, and weakest for defects 
aligned along X-axis. Figure 3b plots the same data in 
the logarithmic coordinates. Rather unexpectedly, one 
can see that over large time intervals the data is consis- 
tent with R oc T 1 ! 2 law, which in turn suggests V oc 
scaling. Up to logarithmic corrections this scaling is what 
one could expect in the static approximation discussed 
above. The same scaling is observed for other PHD con- 
figurations. 

In conclusion, we investigated the motion of penta- 
hepta defects in slightly non-optimal hexagonal patterns 
and their interaction. PHD is stationary only in the per- 
fect pattern with wavenumbers of all three modes equal 
to the onset value (in hexagonal patterns it does not cor- 
respond to the boundary of zig-zag instability). In non- 
optimal hexagonal patterns PHD is driven by the su- 
perposition of two Peach-Kohler forces, corresponding to 
two singular modes. Comparison of the theoretical pre- 
dictions with numerical simulations of (1) showed good 
agreement in both the direction of PHD drift and the 
magnitude of velocity. PHDs attract or repel each other 
depending on the parameter N introduced above. The 
trajectories of interacting defects may be rather compli- 
cated. Furthermore, when two attracting PHDs collide, 
they do not necessarily annihilate, but may give birth to 
another PHD with a different topological structure. 

Equations (1) represent only a simplest possible model 
for hexagonal pattern formation. More realistic models 
derived from first principles (see, e.g., [9]), usually include 
non-variational terms. Nevertheless, we expect that ma- 
jor features of PHD dynamics described in this Letter 
will remain unchanged. We believe that detailed experi- 
ments with thermoconvection or parametric ripples sim- 
ilar to ones which have been performed for dislocations 
in roll patterns [20] could test predictions of our theory. 
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